PredictingThe Onset of the Indian Monsoon

Joint Math Meetings, Seattle, 2025

Sushovan Majhi
Assistant Professor, Data Science
George Washington University, Washington D.C., USA

Today’s Talk

Background on the Indian Monsoon
Topology of a time series
Takens’s Embedding Theorem
Sliding Window
Predicting the onset and withdrawal
Future work

The Indian Monsoon

Rain in India

The climate of India predominately depends on rainfall
Average: 899 mm with variation \pm 20\%
Monsoon is the typical rainy season
- Onset: May-July
- Withdrawal: Sep-Oct
What is the Monsoon?

The Monsoon System

Definition according to the Indian Meteorological Department (IMD)

The seasonal reversal of the direction of winds along the shores of the Indian Ocean, especially in the Arabian Sea, which blow from the southwest for half of the year and from the northeast for the other half.

Only partly understood
The most anticipated weather phenomenon
The biggest financial bet
Notoriously difficult to predict
Affecting the Indian flora and fauna, economy, and agriculture

Monsoon is Coming!

The main objective of the project is develop an early warning system using topological data analysis (TDA).

Climate Science & TDA

Nuisances of Climate Data

high dimensionality
complexity of realistic models
presence of noise
missing values

Successful Applications of TDA

European Topsoil
- Savic, Toth, and Duponchel (2017)
Wildfire
- Kim and Vogel (2019)
Atmospheric River Patters
- Muszynski et al. (2019)
Weather Regimes
- Strommen et al. (2023)

Topology of a Time-series

Chaos

A chaotic system is deterministic but sensitive to initial conditions and manifest unpredictable patterns.

Lorenz System

\begin{aligned} \dot{x} &= \sigma(y-x)\\ \dot{y} &= x(\rho-z)-y\\ \dot{z} &= xy-\beta z \end{aligned} Parameters: \sigma = 10, \rho = 28, \beta = 8 / 3.

Chaos is still discernible when projected on the xy-plane.
How much of it still detectable using only a single signal, e.g. x?

Plot.plot({
  style: { fontFamily: 'monospace', fontSize: 18 },
  marginTop: 50,
  marginBottom: 50,
  height: 2000,
  width: width,
  x: { domain: [-20, 20], ticks: 15 },
  y: { domain: [-23, 30], ticks: 15 },
  marks: [
    Plot.line(lorenzData.slice(0, n), {
      x: "x",
      y: "y",
      strokeWidth: 2,
      stroke: "red"
    }),
    Plot.ruleX([-20]),
    Plot.ruleY([-23])
  ]
});

mutable play = false;
viewof playPause = Inputs.button(" ⏯️", {
  value: false,
  reduce: () => { if(mutable n >= 3999) mutable n = 0; mutable play = !play; } 
});

Time-Delay Embedding

Let’s now consider only one signal: x

Plot.plot({
  style: { fontFamily: 'monospace', fontSize: 18 },
  height: 250,
  width: width,
  x: {label: "t →"},
  y: { ticks: 4 },
  marks: [
    Plot.line(lorenzData.slice(500,600), {
      x: "n",
      y: "x"
    }),
    Plot.ruleX([w, w+delay], { stroke: "red" }),
    Plot.text([[w+0.2, -5], [w+delay+0.6, -5]], {text: [`d₁`, `d₂`]}),
    Plot.ruleX([500]),   
    Plot.ruleY([-12]),   
  ]
});

viewof w = Inputs.range([500, 600-delay-1], { value: 510, step: 1 });

The topology of the time-series is shown below:

Plot.plot({
  style: { fontFamily: 'monospace', fontSize: 18 },
  width: width,
  height: 650,
  x: {label: 'd₁', ticks: 15 },
  y: {label: 'd₂', ticks: 10 },
  marks: [
    Plot.ruleX([22.5]),
    Plot.ruleY([22]),
    Plot.dot(
      delayTS(
        lorenzData.slice(500,600).map((d) => d.z),
        delay
      ),
      {
        stroke: "white",
        r: (d, id) => id == w-500 ? 100 : 20,
        fill: (d, id) => id == w-500 ? "red" : "none"
      }
    )
  ]
});

viewof delay = Inputs.range([1, 50], { value: 1, label: tex`\text{Delay }(\tau)=`, step: 1 });

Sliding Window

We still consider only one signal: x

Plot.plot({
  style: { fontFamily: 'monospace', fontSize: 18 },
  height: 300,
  width: width,
  x: { ticks: 25, label: 't →' },
  y: {ticks: 7 },
  marks: [
    Plot.ruleX([500, 1300]),
    Plot.ruleY([-20]),
    Plot.line(lorenzData.slice(500), {
      x: "n",
      y: "x",
      stroke: (d) => d.n <= 1300 ? "green" : "red"
    })
  ]
});

Below are the time-delay embeddings of the above windows:

Plot.plot({
  style: { fontFamily: 'monospace', fontSize: 22 },
  width: width,
  height: 1500,
  x: { ticks: 15, label: 'd₁'},
  y: { ticks: 15, label: 'd₂'},
  marks: [
    Plot.ruleX([-15]),
    Plot.ruleY([-15]),
    Plot.dot(
      delayTS(
        lorenzData.map((d) => d.x),
        delay1
      ).slice(500, 1300),
      {
        stroke: "green",
        r: 8
      }
    )
  ]
});

Plot.plot({
  style: { fontFamily: 'monospace', fontSize: 22 },
  width: width,
  height: 1500,
  x: { ticks: 15, label: 'd₁' },
  y: { ticks: 15, label: 'd₂'},
  marks: [
    Plot.ruleX([-18]),
    Plot.ruleY([-18]),
    Plot.dot(
      delayTS(
        lorenzData.map((d) => d.x),
        delay1
      ).slice(1300),
      {
        stroke: "red",
        r: 8
      }
    )
  ]
});

viewof delay1 = Inputs.range([1, 50], { value: 1, label: tex`\text{Delay }(\tau)=`, step: 1 });

TDA Pipeline

We use overlapping windows to detect transition to chaos.

take the input signal
set a window size and delay
position the window at the beginning of the time-series, and do the following:
- construct a point-cloud using time-delay embedding
- compute persistence diagram
- additionally, compute summaries like persistence landscapes, L^p norms
slide the window one-step forward, and repeat 3 until the end of the time-series is reached.
In our case, we accumulate the computed norms to output a time-series.

Predicting the Monsoon

The Monsoon Index

Plot.plot({
  style: { fontFamily: 'monospace', fontSize: 18 },
  width: width,
  height: 500,
  marginBottom: 50,
  y: { ticks: 10 },
  marks: [
    Plot.ruleY([-14]),
    Plot.line(
      data,
      {
        x: "Date",
        y: "rain",
      }
    )
  ]
});

data = FileAttachment("idx_1972.csv").csv({typed: true});
norms = FileAttachment("idx_1972_norms.csv").csv({typed: true});

Q: Can you guess the official onset from the signal?
A: June 18 was the official date.

Inference using Sliding Window

Window size =30 and embedding dimension =2.

Plot.plot({
  style: { fontFamily: 'monospace', fontSize: 12 },
  width: width,
  height: 200,
  marks: [
    Plot.ruleY([-14]),
    Plot.dot(
      data,
      {
        x: "Date",
        y: "rain",
        stroke: (d, index) => index >= w1 && index <= w1 + w_size ? "red" : "white"
      }
    )
  ]
});

w_size = 30
viewof w1 = Inputs.range([0, data.length-w_size-1], { value: 0, step: 1 });

Plot.plot({
  style: { fontFamily: 'monospace', fontSize: 10 },
  width: width,
  height: 400,
  x: { domain: [-14,21], label: 'd₁'},
  y: { domain: [-14,21], label: 'd₂'},
  marks: [
    Plot.ruleX([-14]),
    Plot.ruleY([-14]),
    Plot.dot(
      delayTS(
        data.slice(w1, w1 + w_size).map((d) => d.rain),
        delay3
      ),
      {
        stroke: "red",
        r: 4,
      }
    )
  ]
});

viewof delay3 = Inputs.range([1, 20], { value: 1, label: tex`Delay`, step: 1 });

Persistence Pipeline

Go to embedding dimension: 7
Window size: 30 (delay: \tau=1)

The time-series of L^2 norms of persistence landscapes:

Plot.plot({
  style: { fontFamily: 'monospace', fontSize: 13 },
  width: width,
  height: 500,
  marks: [
    Plot.ruleX([Date.parse('1972-01-01')]),
    Plot.ruleX([Date.parse('1972-12-31')], { stroke: 'red' }),
    Plot.ruleY([0]),
    dualAxisY(data, {y: "rain", anchor: "left", color: "white"}),
    dualAxisY(norms, {y: "norm", anchor: "right", color: "red"}),
    Plot.lineY(
      data,
      Plot.normalizeY("extent", {x: "Date", y: "rain", stroke: "white"})
    ),
    Plot.lineY(
      norms,
      Plot.normalizeY("extent", {x: "Date", y: "norm", stroke: "red"})
    )
  ]
});

Onset and Withdrawal of Monsoon

Onset: first occurrence of norm larger than \mu+c_1*\sigma
Withdrawal: last occurrence of norm larger than \mu+c_2*\sigma

c_1=2, c_2=1:

Plot.plot({
  style: { fontFamily: 'monospace', fontSize: 13 },
  width: width,
  height: 500,
  marks: [
    Plot.text([[Date.parse('1972-06-15'), 1], [Date.parse('1972-09-26'), 1]], { text: ['June 15', 'Oct 26'] }),
    Plot.ruleY([0]),
    Plot.ruleX([Date.parse('1972-01-01')]),
    Plot.ruleX([Date.parse('1972-12-31')], { stroke: 'red'}),
    Plot.ruleX([Date.parse('1972-06-15')], { stroke: 'green', strokeWidth: 3 }),
    Plot.ruleX([Date.parse('1972-09-26')], { stroke: 'green', strokeWidth: 3}),
    dualAxisY(data, {y: "rain", anchor: "left", color: "white"}),
    dualAxisY(norms, {y: "norm", anchor: "right", color: "red"}),
    Plot.lineY(
      data,
      Plot.normalizeY("extent", {x: "Date", y: "rain", stroke: "white"})
    ),
    Plot.lineY(
      norms,
      Plot.normalizeY("extent", {x: "Date", y: "norm", stroke: "red"})
    )
  ]
});

Transition to Chaos

Comparing 2-dimensional persistence diagrams

Future Directions

For a working early-warning system, we need a better resolution of the transition period
Enhance interpretability
Instead of relying on delay coordinates, use more than just one singal
Fine-tune the parameters (e.g. window size, embedding dimension, delay) to minimize false-negatives and false-positives
PCA to project the high-dimensional reconstruction of the phase space to gain insights into the attractors
Extend the study to other tropical monsoon systems.

Any Questions?

References

delayTS = function (ts, delay = 1, dim = 2) {
  return ts.map((d, i) => {
    const out = [];
    for (let j = 0; j < dim; j++) out.push(ts[i + j * delay]);
    return out;
  });
};
lorenz = function (
  [x, y, z] = [0, 0, 0],
  n = 100,
  { sigma = 10, r = 28, b = 8 / 3 } = {}
) {
  let F = function (t, [x, y, z]) {
    return [sigma * (y - x), r * x - x * z - y, x * y - b * z];
  };

  let s = new odex.Solver(3);
  let flow = [];

  for (let i = 0; i < n; i++) {
    [x, y, z] = s.solve(F, 0, [x, y, z], 0.01).y;
    flow.push({ n: i, x: x, y: y, z: z });
  }
  return flow;
};
lorenzData = lorenz([1, 1, 1], 4000);
odex = import("https://cdn.skypack.dev/odex@v2.0.4");
mutable n = 0;
{
  let i = 0;
  while (play && i < 4000) {
    yield Promises.delay(3, mutable n++);
  }
};

function dualAxisY(data, {y, ticks = 10, tickFormat, ...options} = {}) {
  const [y1, y2] = d3.extent(Plot.valueof(data, y));
  const scale = d3.scaleLinear().domain([y1, y2]);
  return Plot.axisY(d3.ticks(y1, y2, ticks), {...options, y: scale, tickFormat: scale.tickFormat(ticks, tickFormat)});
}

Kim, Hannah, and Christian Vogel. 2019. “Deciphering Active Wildfires in the Southwestern USA Using Topological Data Analysis.” Climate 7 (12): 135. https://doi.org/10.3390/cli7120135.

Muszynski, Grzegorz, Karthik Kashinath, Vitaliy Kurlin, Michael Wehner, and Prabhat. 2019. “Topological Data Analysis and Machine Learning for Recognizing Atmospheric River Patterns in Large Climate Datasets.” Geoscientific Model Development 12 (2): 613–28. https://doi.org/10.5194/gmd-12-613-2019.

Savic, Aleksandar, Gergely Toth, and Ludovic Duponchel. 2017. “Topological Data Analysis (TDA) Applied to Reveal Pedogenetic Principles of European Topsoil System.” Science of The Total Environment 586 (May): 1091–1100. https://doi.org/10.1016/j.scitotenv.2017.02.095.

Strommen, Kristian, Matthew Chantry, Joshua Dorrington, and Nina Otter. 2023. “A Topological Perspective on Weather Regimes.” Climate Dynamics 60 (5-6): 1415–45. https://doi.org/10.1007/s00382-022-06395-x.